General mixed linear models
Tracy Sánchez
Última actualización: 25 April 2023
Loading Data
[1] 365 3572 7264 8469 6642 5741 4176 7842 9472 9502 7093 4623 2258 5238 1754
[16] 8695 5531 3976 7412 4796 3246 9601 8007 479 8673 2361 3310 9586 2693 4917
[31] 4597Recoding for effects and formating factors
HumanAf$ContextEffectf <-dplyr::recode(HumanAf$Context,
False = -0.5, True= 0.5,
.default = NaN)
HumanAf$AgentPresence <-dplyr::recode(HumanAf$AvatarPresenceCategory,
Omitted = -0.5, Present= 0.5,
.default = NaN)
HumanAf$ContextEffectf <-factor(HumanAf$ContextEffectf,levels= c(-0.5, 0.5),
labels=c('Residential', 'Public'))
HumanAf$AgentPresencef <-factor(HumanAf$AgentPresence,
levels= c(-0.5, 0.5),
labels=c('Omitted', 'Displayed'))Data descriptives
MainVariables <- subset(HumanAf, select = c(AbsolutError, RT))
summary(MainVariables) AbsolutError RT
Min. : 0.00986 Min. : 0.9558
1st Qu.: 12.79062 1st Qu.: 3.7612
Median : 34.99497 Median : 6.3956
Mean : 49.41378 Mean : 8.1636
3rd Qu.: 73.85184 3rd Qu.:10.9006
Max. :179.98339 Max. :29.5264
NA's :1 NA's :1 summary(HumanAf$AgentPresence) Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
-0.5000000 -0.5000000 0.5000000 0.0002036 0.5000000 0.5000000 1 df = HumanAf[complete.cases(HumanAf),]
ggplot(df, aes(x=ContextEffectf, y=AbsolutError, fill=AgentPresencef)) +
geom_boxplot(notch=TRUE,
notchwidth = 0.8,
outlier.colour="red",
outlier.fill="red",
outlier.size=0.5)
ggplot(df, aes(x=ContextEffectf, y=RT, fill=AgentPresencef)) +
geom_boxplot(notch=TRUE,
notchwidth = 0.8,
outlier.colour="red",
outlier.fill="red",
outlier.size=0.5)
library(dplyr)
TwoFactorTable <- HumanAf %>%
group_by(ContextEffectf, AgentPresencef)%>%
summarise(AccuracyMean = mean(AbsolutError, na.rm = TRUE),
n=n(),
AccuracyStandardDev = sd(AbsolutError, na.rm = TRUE),
RTMean = mean(RT, na.rm = TRUE),
RTStandardDev = sd(RT, na.rm = TRUE))`summarise()` has grouped output by 'ContextEffectf'. You can override using
the `.groups` argument.library(tidyr)
Attaching package: 'tidyr'
The following objects are masked from 'package:Matrix':
expand, pack, unpack
The following object is masked from 'package:dlookr':
extractTwoFactorTableUnite <- TwoFactorTable %>%
unite("TwoFactor", ContextEffectf:AgentPresencef, sep= " ", remove = F)
TwoFactorTableUnite <- TwoFactorTableUnite %>%
mutate( AccuracyStandardError=AccuracyStandardDev/sqrt(n)) %>%
mutate( AccuracyStandardIC=AccuracyStandardDev * qt((1-0.05)/2 + .5, n-1)) %>%
mutate( RTStandardError=RTStandardDev/sqrt(n)) Warning in qt((1 - 0.05)/2 + 0.5, n - 1): NaNs producedWarning: Ignoring unknown aesthetics: linetype
Warning: Ignoring unknown aesthetics: linetype
n_distinct(TwoFactorTableUnite$ID)Warning: Unknown or uninitialised column: `ID`.[1] 0Checking for the distribution of Absolut error
df = HumanAf[complete.cases(HumanAf),]
df$AbsolutErrorR <- round(df$AbsolutError, digits = 3)
qqp(df$AbsolutErrorR, "norm")
[1] 2188 594n_distinct(df$ID)[1] 21Linear mixed models
Assessing the need for a multilevel model
interceptOnly <-gls(AbsolutError ~ 1, data = df,
method = "ML")
IDrandomInterceptOnly <-lme(AbsolutError ~ 1, data = df,
random =~1|ID,
method = "ML")
StartlocationsrandomIntercept <-update(IDrandomInterceptOnly, .~.,
random=~1|ID/PointingTaskStartingLocations,
method= "ML")Including Id and starting position as random effects significantly improves the fit of the model
Absolut Error Models
I am adding one main factor at a time
MeaningfulContext <-update(StartlocationsrandomIntercept, .~. + ContextEffectf)
summary(MeaningfulContext)Linear mixed-effects model fit by maximum likelihood
Data: df
AIC BIC logLik
72733.39 72767.69 -36361.7
Random effects:
Formula: ~1 | ID
(Intercept)
StdDev: 15.27947
Formula: ~1 | PointingTaskStartingLocations %in% ID
(Intercept) Residual
StdDev: 12.7486 40.88008
Fixed effects: AbsolutError ~ ContextEffectf
Value Std.Error DF t-value p-value
(Intercept) 50.28221 3.447704 6449 14.584260 0.0000
ContextEffectfPublic -2.23042 1.006543 6449 -2.215926 0.0267
Correlation:
(Intr)
ContextEffectfPublic -0.146
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.1200173 -0.6709893 -0.2286124 0.4965625 3.5612763
Number of Observations: 7038
Number of Groups:
ID PointingTaskStartingLocations %in% ID
21 588 Anova(MeaningfulContext)Analysis of Deviance Table (Type II tests)
Response: AbsolutError
Chisq Df Pr(>Chisq)
ContextEffectf 4.9117 1 0.02668 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Presence <-update(MeaningfulContext, .~. + AgentPresencef)
summary(Presence)Linear mixed-effects model fit by maximum likelihood
Data: df
AIC BIC logLik
72732.14 72773.29 -36360.07
Random effects:
Formula: ~1 | ID
(Intercept)
StdDev: 15.2803
Formula: ~1 | PointingTaskStartingLocations %in% ID
(Intercept) Residual
StdDev: 12.74539 40.87063
Fixed effects: AbsolutError ~ ContextEffectf + AgentPresencef
Value Std.Error DF t-value p-value
(Intercept) 49.38067 3.484115 6448 14.173091 0.0000
ContextEffectfPublic -2.24102 1.006399 6448 -2.226777 0.0260
AgentPresencefDisplayed 1.81276 1.005084 6448 1.803586 0.0713
Correlation:
(Intr) CntxEP
ContextEffectfPublic -0.143
AgentPresencefDisplayed -0.143 -0.006
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.1389159 -0.6672874 -0.2255783 0.4910790 3.5844902
Number of Observations: 7038
Number of Groups:
ID PointingTaskStartingLocations %in% ID
21 588 Anova(Presence)Analysis of Deviance Table (Type II tests)
Response: AbsolutError
Chisq Df Pr(>Chisq)
ContextEffectf 4.9606 1 0.02593 *
AgentPresencef 3.2543 1 0.07124 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1TwofactorInteraction <-update(Presence, .~. + ContextEffectf*AgentPresencef)
summary(TwofactorInteraction)Linear mixed-effects model fit by maximum likelihood
Data: df
AIC BIC logLik
72733.96 72781.97 -36359.98
Random effects:
Formula: ~1 | ID
(Intercept)
StdDev: 15.27345
Formula: ~1 | PointingTaskStartingLocations %in% ID
(Intercept) Residual
StdDev: 12.75142 40.86927
Fixed effects: AbsolutError ~ ContextEffectf + AgentPresencef + ContextEffectf:AgentPresencef
Value Std.Error DF t-value
(Intercept) 49.17619 3.516474 6447 13.984515
ContextEffectfPublic -1.82150 1.414246 6447 -1.287963
AgentPresencefDisplayed 2.22991 1.409627 6447 1.581916
ContextEffectfPublic:AgentPresencefDisplayed -0.84531 2.002603 6447 -0.422108
p-value
(Intercept) 0.0000
ContextEffectfPublic 0.1978
AgentPresencefDisplayed 0.1137
ContextEffectfPublic:AgentPresencefDisplayed 0.6730
Correlation:
(Intr) CntxEP AgntPD
ContextEffectfPublic -0.198
AgentPresencefDisplayed -0.198 0.490
ContextEffectfPublic:AgentPresencefDisplayed 0.138 -0.703 -0.701
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.1449871 -0.6674416 -0.2284260 0.4895499 3.5793448
Number of Observations: 7038
Number of Groups:
ID PointingTaskStartingLocations %in% ID
21 588 Anova(TwofactorInteraction)Analysis of Deviance Table (Type II tests)
Response: AbsolutError
Chisq Df Pr(>Chisq)
ContextEffectf 4.9602 1 0.02594 *
AgentPresencef 3.2544 1 0.07123 .
ContextEffectf:AgentPresencef 0.1783 1 0.67286
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Pairwise<- emmeans(TwofactorInteraction, pairwise ~ ContextEffectf*AgentPresencef)
Pairwise$emmeans
ContextEffectf AgentPresencef emmean SE df lower.CL upper.CL
Residential Omitted 49.2 3.52 20 41.8 56.5
Public Omitted 47.4 3.52 20 40.0 54.7
Residential Displayed 51.4 3.52 20 44.1 58.7
Public Displayed 48.7 3.52 20 41.4 56.1
Degrees-of-freedom method: containment
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
Residential Omitted - Public Omitted 1.821 1.41 6447 1.288 0.5706
Residential Omitted - Residential Displayed -2.230 1.41 6447 -1.582 0.3890
Residential Omitted - Public Displayed 0.437 1.42 6447 0.308 0.9899
Public Omitted - Residential Displayed -4.051 1.43 6447 -2.840 0.0234
Public Omitted - Public Displayed -1.385 1.43 6447 -0.970 0.7668
Residential Displayed - Public Displayed 2.667 1.43 6447 1.871 0.2406
Degrees-of-freedom method: containment
P value adjustment: tukey method for comparing a family of 4 estimates plot(Pairwise[[2]], CIs = TRUE)
library(multcomp);library(multcompView)Loading required package: mvtnormLoading required package: survivalLoading required package: TH.data
Attaching package: 'TH.data'The following object is masked from 'package:MASS':
geyserCLD <- cld(Pairwise,
alpha=0.05,
Letters=letters,
adjust="sidak")I bet you wanted to call this with just object[[1]] - use '[[]]' or which' if I'm wrong.
See '? emm_list' for more informationggplot(CLD,
aes(x = ContextEffectf,
y = emmean,
group = AgentPresencef,
colours = .group)) +
geom_point(aes(shape=AgentPresencef, linetype =AgentPresencef), position=position_dodge(0.3)) +
geom_errorbar(aes(linetype=AgentPresencef,
ymin = lower.CL,
ymax = upper.CL),
position=position_dodge(0.3),
width = 0.2,
size = 0.7) +
theme_bw() +
theme(axis.title = element_text(face = "bold"),
axis.text = element_text(face = "bold"),
plot.caption = element_text(hjust = 0)) +
ylab("Estimated marginal mean\ Absolute angular error") +
xlab("Location location") +
ggtitle ("Marginal Means",
subtitle = "location * Presence") +
labs(caption = paste0(
"Boxes indicate the EM mean. \n",
"Error bars indicate the 95% ",
"confidence interval of the EM mean. \n"),
hjust=0.5) Warning: Ignoring unknown aesthetics: linetype
anova(interceptOnly, IDrandomInterceptOnly, StartlocationsrandomIntercept,
MeaningfulContext, Presence, TwofactorInteraction ) Model df AIC BIC logLik Test
interceptOnly 1 2 73708.98 73722.70 -36852.49
IDrandomInterceptOnly 2 3 72927.44 72948.02 -36460.72 1 vs 2
StartlocationsrandomIntercept 3 4 72736.30 72763.74 -36364.15 2 vs 3
MeaningfulContext 4 5 72733.39 72767.69 -36361.70 3 vs 4
Presence 5 6 72732.14 72773.29 -36360.07 4 vs 5
TwofactorInteraction 6 7 72733.96 72781.97 -36359.98 5 vs 6
L.Ratio p-value
interceptOnly
IDrandomInterceptOnly 783.5359 <.0001
StartlocationsrandomIntercept 193.1437 <.0001
MeaningfulContext 4.9091 0.0267
Presence 3.2536 0.0713
TwofactorInteraction 0.1782 0.6729plot(TwofactorInteraction, which = 1)
plot(MeaningfulContext, which = 1)
GHQ <- glmer(AbsolutError ~ ContextEffectf*AgentPresencef + (1|ID), data = HumanAf,family=gaussian(link = "log"), nAGQ = 25) Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
- Rescale variables?summary(GHQ)Generalized linear mixed model fit by maximum likelihood (Adaptive
Gauss-Hermite Quadrature, nAGQ = 25) [glmerMod]
Family: gaussian ( log )
Formula: AbsolutError ~ ContextEffectf * AgentPresencef + (1 | ID)
Data: HumanAf
AIC BIC logLik deviance df.resid
13423512 13423554 -6711750 13423500 7363
Scaled residuals:
Min 1Q Median 3Q Max
-1.6909 -0.7089 -0.2536 0.4937 3.6907
Random effects:
Groups Name Variance Std.Dev.
ID (Intercept) 264.1 16.25
Residual 1821.6 42.68
Number of obs: 7369, groups: ID, 22
Fixed effects:
Estimate Std. Error t value
(Intercept) 3.8424809 0.0811450 47.353
ContextEffectfPublic -0.0464021 0.0006506 -71.321
AgentPresencefDisplayed 0.0364122 0.0006232 58.430
ContextEffectfPublic:AgentPresencefDisplayed 0.0002670 0.0009040 0.295
Pr(>|z|)
(Intercept) <2e-16 ***
ContextEffectfPublic <2e-16 ***
AgentPresencefDisplayed <2e-16 ***
ContextEffectfPublic:AgentPresencefDisplayed 0.768
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) CntxEP AgntPD
CntxtEffctP -0.004
AgntPrsncfD -0.004 0.499
CntxtEP:APD 0.003 -0.720 -0.690
optimizer (Nelder_Mead) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
- Rescale variables?Anova(GHQ)Analysis of Deviance Table (Type II Wald chisquare tests)
Response: AbsolutError
Chisq Df Pr(>Chisq)
ContextEffectf 10502.1071 1 <2e-16 ***
AgentPresencef 6567.3270 1 <2e-16 ***
ContextEffectf:AgentPresencef 0.0872 1 0.7677
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova(GHQ)Analysis of Variance Table
npar Sum Sq Mean Sq F value
ContextEffectf 1 10613.3 10613.3 5.8264
AgentPresencef 1 6579.6 6579.6 3.6120
ContextEffectf:AgentPresencef 1 0.1 0.1 0.0000emmeans(GHQ, pairwise ~ ContextEffectf:AgentPresencef, type = "response")$emmeans
ContextEffectf AgentPresencef response SE df asymp.LCL asymp.UCL
Residential Omitted 46.6 3.78 Inf 39.8 54.7
Public Omitted 44.5 3.61 Inf 38.0 52.2
Residential Displayed 48.4 3.93 Inf 41.3 56.7
Public Displayed 46.2 3.75 Inf 39.4 54.2
Confidence level used: 0.95
Intervals are back-transformed from the log scale
$contrasts
contrast ratio SE df null z.ratio
Residential Omitted / Public Omitted 1.0475 0.0006815 Inf 1 71.321
Residential Omitted / Residential Displayed 0.9642 0.0006009 Inf 1 -58.430
Residential Omitted / Public Displayed 1.0098 0.0006456 Inf 1 15.207
Public Omitted / Residential Displayed 0.9205 0.0005872 Inf 1 -129.829
Public Omitted / Public Displayed 0.9640 0.0006305 Inf 1 -56.077
Residential Displayed / Public Displayed 1.0472 0.0006569 Inf 1 73.548
p.value
<.0001
<.0001
<.0001
<.0001
<.0001
<.0001
P value adjustment: tukey method for comparing a family of 4 estimates
Tests are performed on the log scale plot(fitted(GHQ), residuals(GHQ), xlab = "Fitted Values", ylab = "Residuals")
abline(h = 0, lty = 2)
lines(smooth.spline(fitted(GHQ), residuals(GHQ)))
interceptOnly <-gls(log(AbsolutError) ~ 1, data = df,
method = "ML")
IDrandomInterceptOnly <-lme(log(AbsolutError) ~ 1, data = df,
random =~1|ID,
method = "ML")
StartlocationsrandomIntercept <-update(IDrandomInterceptOnly, .~.,
random=~1|ID/PointingTaskStartingLocations,
method= "ML")
MeaningfulContext <-update(StartlocationsrandomIntercept, .~. + ContextEffectf)
Presence <-update(MeaningfulContext, .~. + AgentPresencef)
TwofactorInteraction <-update(Presence, .~. + ContextEffectf*AgentPresencef)
summary(TwofactorInteraction)Linear mixed-effects model fit by maximum likelihood
Data: df
AIC BIC logLik
23851.32 23899.34 -11918.66
Random effects:
Formula: ~1 | ID
(Intercept)
StdDev: 0.4755101
Formula: ~1 | PointingTaskStartingLocations %in% ID
(Intercept) Residual
StdDev: 0.306225 1.281173
Fixed effects: log(AbsolutError) ~ ContextEffectf + AgentPresencef + ContextEffectf:AgentPresencef
Value Std.Error DF
(Intercept) 3.262988 0.10900052 6447
ContextEffectfPublic -0.069948 0.04406511 6447
AgentPresencefDisplayed 0.089385 0.04394283 6447
ContextEffectfPublic:AgentPresencefDisplayed -0.048219 0.06238529 6447
t-value p-value
(Intercept) 29.935528 0.0000
ContextEffectfPublic -1.587387 0.1125
AgentPresencefDisplayed 2.034123 0.0420
ContextEffectfPublic:AgentPresencefDisplayed -0.772919 0.4396
Correlation:
(Intr) CntxEP AgntPD
ContextEffectfPublic -0.199
AgentPresencefDisplayed -0.199 0.491
ContextEffectfPublic:AgentPresencefDisplayed 0.139 -0.704 -0.703
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-5.8868056 -0.4309144 0.1913041 0.6743222 2.2573072
Number of Observations: 7038
Number of Groups:
ID PointingTaskStartingLocations %in% ID
21 588 Anova(TwofactorInteraction)Analysis of Deviance Table (Type II tests)
Response: log(AbsolutError)
Chisq Df Pr(>Chisq)
ContextEffectf 9.0114 1 0.002683 **
AgentPresencef 4.3937 1 0.036071 *
ContextEffectf:AgentPresencef 0.5977 1 0.439440
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Pairwise<- emmeans(TwofactorInteraction, pairwise ~ ContextEffectf*AgentPresencef, type='response')
Pairwise $emmeans
ContextEffectf AgentPresencef response SE df lower.CL upper.CL
Residential Omitted 26.1 2.85 20 20.8 32.8
Public Omitted 24.4 2.66 20 19.4 30.6
Residential Displayed 28.6 3.12 20 22.8 35.9
Public Displayed 25.4 2.77 20 20.2 31.9
Degrees-of-freedom method: containment
Confidence level used: 0.95
Intervals are back-transformed from the log scale
$contrasts
contrast ratio SE df null t.ratio
Residential Omitted / Public Omitted 1.072 0.0473 6447 1 1.587
Residential Omitted / Residential Displayed 0.914 0.0402 6447 1 -2.034
Residential Omitted / Public Displayed 1.029 0.0454 6447 1 0.653
Public Omitted / Residential Displayed 0.853 0.0379 6447 1 -3.588
Public Omitted / Public Displayed 0.960 0.0426 6447 1 -0.927
Residential Displayed / Public Displayed 1.125 0.0499 6447 1 2.667
p.value
0.3859
0.1755
0.9145
0.0019
0.7902
0.0384
Degrees-of-freedom method: containment
P value adjustment: tukey method for comparing a family of 4 estimates
Tests are performed on the log scale plot(Pairwise[[2]])
anova(interceptOnly, IDrandomInterceptOnly, StartlocationsrandomIntercept,
MeaningfulContext, Presence, TwofactorInteraction ) Model df AIC BIC logLik Test
interceptOnly 1 2 24731.54 24745.26 -12363.77
IDrandomInterceptOnly 2 3 23939.25 23959.83 -11966.62 1 vs 2
StartlocationsrandomIntercept 3 4 23859.22 23886.66 -11925.61 2 vs 3
MeaningfulContext 4 5 23852.31 23886.61 -11921.16 3 vs 4
Presence 5 6 23849.92 23891.07 -11918.96 4 vs 5
TwofactorInteraction 6 7 23851.32 23899.34 -11918.66 5 vs 6
L.Ratio p-value
interceptOnly
IDrandomInterceptOnly 794.2943 <.0001
StartlocationsrandomIntercept 82.0308 <.0001
MeaningfulContext 8.9072 0.0028
Presence 4.3924 0.0361
TwofactorInteraction 0.5971 0.4397plot(TwofactorInteraction, which = 1)
model <- glmmPQL(RT ~ ContextEffectf*AgentPresencef, ~1|ID/PointingTaskStartingLocations, family = gaussian(link = "log"), data = HumanAf, verbose = FALSE)
summary(model)Linear mixed-effects model fit by maximum likelihood
Data: HumanAf
AIC BIC logLik
NA NA NA
Random effects:
Formula: ~1 | ID
(Intercept)
StdDev: 0.2478035
Formula: ~1 | PointingTaskStartingLocations %in% ID
(Intercept) Residual
StdDev: 0.2197559 5.094753
Variance function:
Structure: fixed weights
Formula: ~invwt
Fixed effects: RT ~ ContextEffectf * AgentPresencef
Value Std.Error DF
(Intercept) 2.0639381 0.05550721 6750
ContextEffectfPublic 0.0276047 0.01977188 6750
AgentPresencefDisplayed -0.0097499 0.01986404 6750
ContextEffectfPublic:AgentPresencefDisplayed -0.0035024 0.02808134 6750
t-value p-value
(Intercept) 37.18324 0.0000
ContextEffectfPublic 1.39616 0.1627
AgentPresencefDisplayed -0.49083 0.6236
ContextEffectfPublic:AgentPresencefDisplayed -0.12472 0.9007
Correlation:
(Intr) CntxEP AgntPD
ContextEffectfPublic -0.178
AgentPresencefDisplayed -0.176 0.490
ContextEffectfPublic:AgentPresencefDisplayed 0.124 -0.702 -0.703
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.1014455 -0.6668605 -0.3148692 0.3570815 4.4095516
Number of Observations: 7369
Number of Groups:
ID PointingTaskStartingLocations %in% ID
22 616 Response Time models
Checking for the distribution of RT
df = HumanAf[complete.cases(HumanAf),]
df$RTr <- round(df$RT, digits = 3)
qqp(df$RT, "norm")
[1] 6521 6473qqp(df$RT, "lnorm")
[1] 6521 6473interceptOnlyt <-gls(log(RTr) ~ 1, data = df,
method = "ML")
IDrandomInterceptOnlyt <-lme(log(RTr) ~ 1, data = df,
random =~1|ID,
method = "ML")
StartlocationsrandomInterceptt <-lme(log(RTr) ~ 1, data = df,
random=~1|ID|PointingTaskStartingLocations,
method= "ML")
MeaningfulContext <-update(StartlocationsrandomInterceptt, .~. + ContextEffectf)
Presence <-update(MeaningfulContext, .~. + AgentPresencef)
TwofactorInteraction <-update(Presence, .~. + ContextEffectf*AgentPresencef)
summary(TwofactorInteraction)Linear mixed-effects model fit by maximum likelihood
Data: df
AIC BIC logLik
14903.06 14957.94 -7443.532
Random effects:
Formula: ~1 | ID | PointingTaskStartingLocations
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.09028261 (Intr)
1 | IDTRUE 0.09028261 -0.676
Residual 0.69492859
Fixed effects: log(RTr) ~ ContextEffectf + AgentPresencef + ContextEffectf:AgentPresencef
Value Std.Error DF
(Intercept) 1.8835518 0.02146533 7007
ContextEffectfPublic -0.0008648 0.02351088 7007
AgentPresencefDisplayed -0.0423609 0.02345250 7007
ContextEffectfPublic:AgentPresencefDisplayed 0.0210118 0.03320587 7007
t-value p-value
(Intercept) 87.74857 0.0000
ContextEffectfPublic -0.03678 0.9707
AgentPresencefDisplayed -1.80624 0.0709
ContextEffectfPublic:AgentPresencefDisplayed 0.63277 0.5269
Correlation:
(Intr) CntxEP AgntPD
ContextEffectfPublic -0.540
AgentPresencefDisplayed -0.540 0.494
ContextEffectfPublic:AgentPresencefDisplayed 0.381 -0.707 -0.706
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-2.659401131 -0.736408210 -0.008763436 0.746873336 2.287183588
Number of Observations: 7038
Number of Groups: 28 Anova(TwofactorInteraction)Analysis of Deviance Table (Type II tests)
Response: log(RTr)
Chisq Df Pr(>Chisq)
ContextEffectf 0.3366 1 0.56177
AgentPresencef 3.6869 1 0.05484 .
ContextEffectf:AgentPresencef 0.4006 1 0.52676
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Pairwiset<- emmeans(TwofactorInteraction, pairwise ~ ContextEffectf*AgentPresencef, type='response')
Pairwiset $emmeans
ContextEffectf AgentPresencef response SE df lower.CL upper.CL
Residential Omitted 6.58 0.141 27 6.29 6.87
Public Omitted 6.57 0.142 27 6.29 6.87
Residential Displayed 6.30 0.136 27 6.03 6.59
Public Displayed 6.43 0.138 27 6.15 6.72
Degrees-of-freedom method: containment
Confidence level used: 0.95
Intervals are back-transformed from the log scale
$contrasts
contrast ratio SE df null t.ratio
Residential Omitted / Public Omitted 1.00 0.0235 7007 1 0.037
Residential Omitted / Residential Displayed 1.04 0.0245 7007 1 1.806
Residential Omitted / Public Displayed 1.02 0.0239 7007 1 0.950
Public Omitted / Residential Displayed 1.04 0.0246 7007 1 1.756
Public Omitted / Public Displayed 1.02 0.0240 7007 1 0.908
Residential Displayed / Public Displayed 0.98 0.0230 7007 1 -0.858
p.value
1.0000
0.2704
0.7777
0.2950
0.8007
0.8267
Degrees-of-freedom method: containment
P value adjustment: tukey method for comparing a family of 4 estimates
Tests are performed on the log scale ref_grid(TwofactorInteraction)'emmGrid' object with variables:
ContextEffectf = Residential, Public
AgentPresencef = Omitted, Displayed
Transformation: "log" plot(Pairwiset[[2]])
anova(interceptOnlyt, IDrandomInterceptOnlyt, StartlocationsrandomInterceptt,
MeaningfulContext, Presence, TwofactorInteraction ) Model df AIC BIC logLik Test
interceptOnlyt 1 2 14934.25 14947.97 -7465.124
IDrandomInterceptOnlyt 2 3 13601.46 13622.04 -6797.730 1 vs 2
StartlocationsrandomInterceptt 3 5 14901.46 14935.76 -7445.731 2 vs 3
MeaningfulContext 4 6 14903.15 14944.30 -7445.575 3 vs 4
Presence 5 7 14901.46 14949.48 -7443.732 4 vs 5
TwofactorInteraction 6 8 14903.06 14957.94 -7443.532 5 vs 6
L.Ratio p-value
interceptOnlyt
IDrandomInterceptOnlyt 1334.7889 <.0001
StartlocationsrandomInterceptt 1296.0027 <.0001
MeaningfulContext 0.3126 0.5761
Presence 3.6852 0.0549
TwofactorInteraction 0.4005 0.5268plot(TwofactorInteraction, which = 1)
CLD <- cld(Pairwiset,
alpha=0.05,
Letters=letters,
adjust="sidak")I bet you wanted to call this with just object[[1]] - use '[[]]' or which' if I'm wrong.
See '? emm_list' for more informationggplot(CLD,
aes(x = ContextEffectf,
y = response,
group = AgentPresencef,
colours = .group)) +
geom_point(aes(shape=AgentPresencef, color=AgentPresencef), position=position_dodge(0.3)) +
geom_errorbar(aes(color=AgentPresencef,
ymin = lower.CL,
ymax = upper.CL),
position=position_dodge(0.3),
width = 0.2,
size = 0.7) +
theme_bw() +
theme(axis.title = element_text(face = "bold"),
axis.text = element_text(face = "bold"),
plot.caption = element_text(hjust = 0)) +
ylab("Estimated marginal mean\ Response Time in Seconds") +
xlab("Location location") +
ggtitle ("Marginal Means",
subtitle = "location * Presence") +
labs(caption = paste0(
"Boxes indicate the EM mean. \n",
"Error bars indicate the 95% ",
"confidence interval of the EM mean. \n"),
hjust=0.5) 
HumanAf$Agent_Category <- with(HumanAf, ave(seq_along(ID), ID, FUN = function(x) sample(c(rep('Action', ceiling(length(x)*0.6)), rep('Standing', length(x) - ceiling(length(x)*0.6))))))library(dplyr)
TwoFactorTable <- HumanAf %>%
group_by(ContextEffectf, AgentPresencef, Agent_Category)%>%
summarise(AccuracyMean = mean(AbsolutError, na.rm = TRUE),
n=n(),
AccuracyStandardDev = sd(AbsolutError, na.rm = TRUE),
RTMean = mean(RT, na.rm = TRUE),
RTStandardDev = sd(RT, na.rm = TRUE))`summarise()` has grouped output by 'ContextEffectf', 'AgentPresencef'. You can
override using the `.groups` argument.library(tidyr)
TwoFactorTableUnite <- TwoFactorTable %>%
unite("TwoFactor", ContextEffectf:Agent_Category, sep= " ", remove = F)
TwoFactorTableUnite <- TwoFactorTableUnite %>%
mutate( AccuracyStandardError=AccuracyStandardDev/sqrt(n)) %>%
mutate( AccuracyStandardIC=AccuracyStandardDev * qt((1-0.05)/2 + .5, n-1)) %>%
mutate( RTStandardError=RTStandardDev/sqrt(n)) Warning in qt((1 - 0.05)/2 + 0.5, n - 1): NaNs producedggplot(data=subset(TwoFactorTableUnite, !is.na(AgentPresencef)),
aes(x = ContextEffectf,
y = AccuracyMean,
group = Agent_Category)) +
geom_point(aes(shape=Agent_Category, linetype =Agent_Category), position=position_dodge(0.3)) +
geom_errorbar(aes(linetype=Agent_Category,
ymin= AccuracyMean-AccuracyStandardError,
ymax=AccuracyMean+AccuracyStandardError),
position=position_dodge(0.3),
width = 0.2,
size = 0.7) +
facet_wrap(~ AgentPresencef) +
theme_bw() +
theme(axis.title = element_text(face = "bold"),
axis.text = element_text(face = "bold"),
plot.caption = element_text(hjust = 0)) +
ylab("Absolute angular error in degrees") +
xlab("Location location") +
ggtitle ("Accuracy performance",
subtitle = "The effect of location and presence") +
labs(caption = paste0(
"Error bars indicate one Standard Error \n"),
hjust=0.5) Warning: Ignoring unknown aesthetics: linetype